Análisis de regresión lineal

Autor

Omar E. Barrantes Sotela

Fecha de Publicación

7 de septiembre de 2022

Introducción

El análisis de regresión es un grupo de procesos estadísticos para determinar la relación entre las variables del conjunto de datos. Generalmente, el análisis de regresión se utiliza para determinar la relación entre las variables explicativas y respuesta del conjunto de datos.

El análisis de regresión ayuda a comprender cómo cambian las variables respuesta cuando una de las variables explicativas cambia y otras variables se mantienen constantes. Esto ayuda a construir un modelo de regresión y, además, ayuda a pronosticar los valores con respecto a un cambio en una de las variables explicativas.

Existen 4 tipos de técnicas de análisis de regresión:

  • regresión lineal,

  • regresión logística,

  • regresión logística multinomial y

  • regresión logística ordinal.

En este caso se concentrará la atención en el primer tipo.

Demostración

Contexto de los datos

Para este ejemplo, se usan los datos de una empresa aseguradora. Consiste en una base de datos de 1338 registros de personas (descarga de los datos).

Formato

Consiste de 7 variables:

Variable Descripción
age Edad del asegurado.
sex Sexo del asegurado.
bmi Índice de Masa Corporal de la persona.
children Cantidad de hijos asegurados por el seguro médico.
smoker Indicador si la persona fuma.
region Región donde vive el asegurado en E.E.U.U.: northeast, southeast, southwest, northwest.
charges Costos médicos individuales facturados por el seguro de salud.

Se procede a cargar los datos.

library(readr)
#url <- "https://drive.google.com/file/d/1EOSYAQqksn70clDAut3FqLbj0a0GMCM4/view?usp=sharing"
#download.file(url, destfile = "insurance.csv")
insurance <- read_csv("insurance.csv")
Rows: 1338 Columns: 7
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (3): sex, smoker, region
dbl (4): age, bmi, children, charges

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

AED

Se realiza una AED para revisar las características de los datos a utilizar en el análisis.

  • Resumen estadístico
summary(insurance)
      age            sex                 bmi           children    
 Min.   :18.00   Length:1338        Min.   :15.96   Min.   :0.000  
 1st Qu.:27.00   Class :character   1st Qu.:26.30   1st Qu.:0.000  
 Median :39.00   Mode  :character   Median :30.40   Median :1.000  
 Mean   :39.21                      Mean   :30.66   Mean   :1.095  
 3rd Qu.:51.00                      3rd Qu.:34.69   3rd Qu.:2.000  
 Max.   :64.00                      Max.   :53.13   Max.   :5.000  
    smoker             region             charges     
 Length:1338        Length:1338        Min.   : 1122  
 Class :character   Class :character   1st Qu.: 4740  
 Mode  :character   Mode  :character   Median : 9382  
                                       Mean   :13270  
                                       3rd Qu.:16640  
                                       Max.   :63770

Se procede a revisar la distribución de las variables edad y IMC.

library(ggplot2)
ggplot(insurance, aes(x=age)) + 
  geom_histogram(binwidth = 10, 
                 colour="blue", 
                 fill="#5588EE") +
  labs(x= "Edad", y="Frecuencia" ) +
  theme_classic()

ggplot(insurance, aes(x=bmi)) + 
  geom_histogram(binwidth = 5, 
                 colour="#E69F00", 
                 fill="#FEBA4F") +
  labs(x= "IMC", y="Frecuencia" ) +
  theme_classic()

Modelo de regresión lineal

Se utilizan las siguientes librerías para el análisis de regresión.

library(car)      # funciones para el ajuste de modelos de regresión
library(visreg)   # visualiza regresiones
library(lmtest)   # prueba modelos lineales
library(corrplot) # gráfico correlaciones
library(corrgram) # más gráficos de correlaciones

Primer modelo se trata de asociar las variables edad e IMC.

library(knitr)
reg01 <- lm(age ~ bmi, data=insurance)
summary(reg01)

Call:
lm(formula = age ~ bmi, data = insurance)

Residuals:
     Min       1Q   Median       3Q      Max 
-26.8632 -12.5477   0.3053  11.6915  26.7248 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 31.48728    1.95884  16.074  < 2e-16 ***
bmi          0.25176    0.06266   4.018 6.19e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 13.97 on 1336 degrees of freedom
Multiple R-squared:  0.01194,   Adjusted R-squared:  0.0112 
F-statistic: 16.15 on 1 and 1336 DF,  p-value: 6.194e-05
kable(reg01$coefficients)
Tabla 1: Coeficentes de la correlación entre la edad y IMC.
x
(Intercept) 31.4872787
bmi 0.2517577

Las estimaciones para los parámetros del modelo: el valor de la intersección \(( b= 31.4872787)\) y el efecto estimado de la edad sobre el IMC \((m=0.2517577)\).

El resultados del coeficiente de correlación \(R=0.1092719\) y el coeficiente de determinación \(R^2 =0.0119403\) sugieren un asociación muy baja o nula de las variables.

Preguntas generadoras

  1. ¿Qué nos indica la baja asociación de las variables edad y imc?

  2. ¿Es posible otro tipo de relación entre otras variables?

Estratificación o niveles

La creación de estratos o niveles puede ayudar a observar patrones de asociación de variables, que de otra manera son difíciles de notar. En este caso se examina entre subconjuntos según sexo y la condición de fumador del asegurado. La asociación de variables es IMC con respecto a el costo de atención médica.

hist.col <- c("#999999", "#E69F00", "#56B4E9")

ggplot(insurance, aes(x=bmi, fill=sex, color=sex)) +
  geom_histogram(binwidth = 5, position="identity", alpha=0.5) +
  labs(x= "IMC", y="Frecuencia" ) +
  scale_color_manual(values=hist.col)+
  scale_fill_manual(values=hist.col) +
  theme_bw()

ggplot(insurance, aes(x=bmi, fill=smoker, color=smoker)) +
  geom_histogram(binwidth = 5,position="identity", alpha=0.5) +
  labs(x= "IMC", y="Frecuencia" ) +
  scale_color_manual(values= hist.col)+
  scale_fill_manual(values= hist.col) +
  theme_bw()

ggplot(insurance, aes(x=charges, fill=sex, color=sex)) +
  geom_histogram(binwidth = 5000,position="identity", alpha=0.5) +
  labs(x= "Prima del seguro médico", y="Frecuencia" ) +
  scale_color_manual(values= hist.col)+
  scale_fill_manual(values= hist.col) +
  theme_bw()

ggplot(insurance, aes(x=charges, fill=smoker, color=smoker)) +
  geom_histogram(binwidth = 5000,position="identity", alpha=0.5) +
  labs(x= "Prima del seguro médico", y="Frecuencia" ) +
  scale_color_manual(values= hist.col)+
  scale_fill_manual(values= hist.col) +
  theme_bw()

ggplot(insurance, aes(x=bmi, y=charges, fill=smoker)) +
  geom_boxplot()+
  scale_fill_manual(values= hist.col) +
  labs(x="IMC", y="Costos médicos") +
  theme_bw()

ggplot(insurance, aes(x=bmi, y=charges, fill=sex)) +
  geom_boxplot()+
  scale_fill_manual(values= hist.col) +
  labs(x="IMC", y="Costos médicos") +
  theme_bw()

Análisis de correlación por estratos

library(dplyr)
library(tidyr)
mat.corr <- insurance |> 
  group_by(sex,smoker) |>
    summarise(r=cor(bmi,charges))
kable(mat.corr)
Tabla 2: Correlación de las variables IMC ~ Costos médicos individuales agrupado por su condición de fumador y el sexo de la persona.
sex smoker r
female no 0.0750925
female yes 0.8459098
male no 0.0959043
male yes 0.7693554

También es posible obtener los coeficientes del modelo lineal según estrato.

library(broom)
reg02 <- insurance |> 
    group_by(sex,smoker) |>
    do(tidy(lm(charges~bmi, .)))
kable(reg02)
Tabla 3: Coeficientes de correlación de las variables IMC ~ costos médicos individuales agrupado por su condición de fumador y el sexo de la persona.
sex smoker term estimate std.error statistic p.value
female no (Intercept) 6407.22480 1364.35924 4.696142 0.0000034
female no bmi 77.11556 43.86509 1.758017 0.0793060
female yes (Intercept) -14075.51629 2720.22707 -5.174390 0.0000010
female yes bmi 1511.55493 89.65063 16.860505 0.0000000
male no (Intercept) 5271.76507 1313.41927 4.013772 0.0000686
male no bmi 91.49778 41.84682 2.186493 0.0292289
male yes (Intercept) -12577.62840 3076.25865 -4.088612 0.0000691
male yes bmi 1448.05010 95.95906 15.090290 0.0000000

Se construye una tabla más compleja con los valores de los parámetros de los modelos creados a partir de la estratificación.

reg02b <- insurance |> 
    group_by(sex,smoker) |>
    do(tidy(cor.test(~charges+bmi, data = . )))

kable(reg02b)
Tabla 4: Correlación de las variables IMC ~ Costos médicos individuales agrupado por su condición de fumador y el sexo de la persona.
sex smoker estimate statistic p.value parameter conf.low conf.high method alternative
female no 0.0750925 1.758017 0.0793060 545 -0.0087984 0.1579338 Pearson’s product-moment correlation two.sided
female yes 0.8459098 16.860505 0.0000000 113 0.7842804 0.8910077 Pearson’s product-moment correlation two.sided
male no 0.0959043 2.186493 0.0292289 515 0.0097493 0.1806458 Pearson’s product-moment correlation two.sided
male yes 0.7693554 15.090290 0.0000000 157 0.6971960 0.8260810 Pearson’s product-moment correlation two.sided

Identificación de las variables más cor-relacionadas y nuevo modelo de regresión

Con los datos anteriores (Tabla 2, Tabla 3 y Tabla 4), es posible confirmar que la condición de fumado en las personas implica una mayor asociación entre el IMC y los costos médicos individuales. Por lo cual es importante detallar más en esta relación.

insurance |> 
  drop_na(bmi,age) |> 
  filter(smoker == "yes") |> 
  cor.test(~charges+bmi, data = _ )

    Pearson's product-moment correlation

data:  charges and bmi
t = 22.496, df = 272, p-value < 2.2e-16
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.7606770 0.8442932
sample estimates:
      cor 
0.8064806 
library(ggplot2)
library(ggpubr)

dat.smoker <- insurance |>
  drop_na(bmi,age) |> 
  filter(smoker == "yes")

reg03.plot <- ggplot(dat.smoker, aes(x=bmi, y=charges))+
                     geom_point() +
  geom_smooth(method="lm", col="blue") +
  stat_regline_equation(label.x = 40, label.y = 20000)+
  theme_bw() +
  labs(x = "IMC", y="Prima del seguro médico en dólares")
reg03.plot
`geom_smooth()` using formula 'y ~ x'

reg03 <- lm(charges ~ bmi, data=dat.smoker)

Las estimaciones para los parámetros del modelo 3: el valor de la intersección \(( b= -1.3186576\times 10^{4})\) y el efecto estimado del IMC sobre el costo médico individual de los asegurados \((m=1473.1062547)\).

El resultados del coeficiente de correlación \(R=0.8064806\) y el coeficiente de determinación \(R^2 =0.650411\) sugieren un asociación muy alta de las variables.

Gráficos de diagnóstico

Es una buena idea producir también gráficos de diagnóstico para analizar los residuos del modelo y asegurarnos de que el modelo lineal sea apropiado para usar para los datos particulares con los que estamos trabajando.

influenceIndexPlot(reg03,id=TRUE)

Gráficos de influencia

influencePlot(reg03, id=list(method="identify"))

Predicción de valores

Que sucedería con mis costos médicos individuales, sí fumo y además puedo presentar los siguientes IMC: 33, 28, 20.

prediccion <- predict(reg03, newdata = data.frame(bmi=c(33,28,20)))
prediccion
       1        2        3 
35425.93 28060.40 16275.55